library(readr)
library(dplyr)
library(tidyr)
library(forcats)
library(stringr)
library(ggplot2)
library(ggtext)
library(lubridate)
library(purrr)
library(sf)
library(terra)
library(rnaturalearth)
library(rcartocolor)
library(colorspace)
library(patchwork)
library(amt)
library(corrr)
library(here)
## plot theme
theme_set(d6::theme_d6(base_size = 18))
theme_update(
panel.grid.major.x = element_blank(),
panel.grid.major.y = element_line(color = "#fefefe", linewidth = .5),
panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "grey95", color = "transparent"),
panel.border = element_rect(fill = "transparent", color = "transparent"),
strip.text = element_markdown(color = "black", hjust = .5, face = "bold", size = rel(1.3)),
legend.text = element_text(family = "PT Mono"),
legend.position = "top",
legend.box = "vertical",
legend.spacing = unit(.2, "lines"),
axis.text = element_text(family = "PT Mono"),
axis.text.x = element_text(margin = margin(t = 5)),
axis.text.y = element_text(margin = margin(r = 5)),
axis.ticks.length = unit(0, "lines"),
panel.spacing.x = unit(.7, "lines"),
panel.spacing.y = unit(.9, "lines"),
plot.margin = margin(rep(.5, 4))
)
## path to processed data
dir_proc <- here("output", "empirical", "data_proc")
## path to analysis results
dir_est <- here("output", "empirical", "results")
## path + folder for movement plots
dir_move <- here("plots", "emp_tracks")
if (!dir.exists(dir_move)) dir.create(dir_move, recursive = TRUE)
## path + folder for plots
dir_plots <- here("plots")
if (!dir.exists(dir_plots)) dir.create(dir_plots, recursive = TRUE)The 4th Dimension of Movement
Analysis (empirical tracks)
Preparation
Data Import
Movement Data
df_boars <- read_rds(paste0(dir_proc, "/boars_resampled_3035.rds"))
sf_boars <-
df_boars %>%
filter(!is.na(x), !is.na(y)) %>%
st_as_sf(
coords = c("x", "y"),
crs = "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs"
) %>%
st_make_valid() %>%
st_transform(crs = 3035) %>%
mutate(
x = st_coordinates(.)[,1],
y = st_coordinates(.)[,2],
t = ymd_hms(time),
id = str_replace(id, "_", " "),
res_label = paste(res, "minutes"),
res_label = fct_reorder(res_label, res)
)
## color palette for boars
n <- length(unique(df_boars$id))
pal <- carto_pal(name = "Bold", n = n+1)[1:n]
names(pal) <- unique(df_boars$id)Landscape Variables
## turn into data frame for plotting
ras_env <- rast(here("output", "empirical", "data_proc", "env_vars_100m.tif"))
df_env <-
ras_env %>%
terra::as.data.frame(xy = TRUE) %>%
as_tibble() %>%
mutate(hab = factor(hab, levels = c(0, 1), labels = c(" non-forest ", " forest")))
bln <- sf::st_transform(d6berlin::sf_berlin, crs = st_crs(ras_env))Figure 2
Movement tracks of the wild boars
sf_boars_move_12h <-
sf_boars |>
filter(res == 720) |>
group_by(id) |>
mutate(t_id = row_number()) |>
arrange(time)
sf_boars_centr_12h <-
sf_boars_move_12h %>%
group_by(id) %>%
summarize(geometry = st_union(geometry)) %>%
st_centroid() %>%
mutate(x = st_coordinates(.)[,1], y = st_coordinates(.)[,2])
set.seed(1)
## Plot
p_1 <-
ggplot(df_env) +
geom_tile(aes(x, y, fill = imp)) +
geom_sf(data = bln, fill = NA, linewidth = 2.7, color = "white") +
geom_sf(data = bln, fill = NA, linewidth = .9, color = "#415c27") +
ggforce::geom_mark_hull(
data = sf_boars_centr_12h,
aes(x, y, group = id, label = paste0("WB", id)),
color = "transparent", con.colour = "grey30",
label.colour = "grey20", label.family = "PT Sans Narrow", label.fontsize = 20,
label.buffer = unit(12, 'mm'), label.margin = margin(2, 2, 1, 2, 'mm'),
expand = unit(1, 'mm'), con.cap = unit(.01, 'mm'), con.size = .7
) +
geom_sf(
data = sf_boars_move_12h,
aes(group = id),
inherit.aes = FALSE,
color = "black",
size = 1.3
) +
geom_path(
data = sf_boars_move_12h,
aes(x, y, group = id),
inherit.aes = FALSE,
linewidth = .8
) +
geom_sf(
data = sf_boars_move_12h,
aes(color = t_id, group = id),
size = .7,
stroke = .3
) +
scale_color_carto_c(palette = "PinkYl", name = "Movement step (relative)") +
scale_fill_gradient(low = "grey95", high = "grey40", name = "Imperviousness (scaled)") +
coord_sf(expand = FALSE, clip = "off") +
guides(
fill = guide_colorbar(
order = 2,
barwidth = unit(24, "lines"), barheight = unit(.6, "lines"),
title.position = "top", title.hjust = .5
),
color = guide_colorbar(
order = 1,
barwidth = unit(24, "lines"), barheight = unit(.6, "lines"),
title.position = "top", title.hjust = .5)
) +
labs(x = "Longitude", y = "Latitude") +
theme(
plot.title = element_text(hjust = .5),
legend.box = "horizontal",
legend.box.margin = margin(t = -5)
)
ggsave(paste0(dir_plots, "/2_emp_boars_all_12h.png"), width = 12, height = 8, dpi = 600, bg = "white")
p_1Single tracks per ID and sampling interval
plot_move_map <- function(boar) {
sf_id <-
sf_boars |>
filter(id == boar) |>
add_count(res) |>
mutate(label = paste0(res, " minutes<br><span style='font-weight:400;font-size:18pt;'>(", n, " locations)</span>")) |>
arrange(res, time) |>
mutate(label = fct_inorder(label))
bbox <-
sf_id|>
st_buffer(1000) |>
st_bbox()
## Plot
g <-
ggplot(df_env, aes(x, y)) +
geom_tile(aes(fill = imp)) +
geom_sf(data = sf_id, color = "grey10", size = 1.7, stroke = .3) +
geom_path(data = sf_id, aes(x, y, group = boar), inherit.aes = FALSE,
color = "grey10", linewidth = .3) +
geom_sf(data = sf_id, aes(color = boar), size = .7) +
coord_sf(xlim = c(bbox$xmin, bbox$xmax),
ylim = c(bbox$ymin, bbox$ymax), expand = FALSE) +
facet_wrap(~label, nrow = 1) +
scale_color_manual(values = pal, guide = "none") +
scale_fill_gradient(low = "grey95", high = "grey40", guide = "none") +
labs(x = NULL, y = NULL, title = paste("Wild boar", boar)) +
theme(
plot.title = element_text(size = rel(2)),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
panel.spacing.x = unit(2.3, "lines")
)
ggsave(paste0(dir_move, "/trajectory_boar_", boar, ".png"),
width = 25, height = 8.5, dpi = 600, bg = "white")
return(g)
}
if(!file.exists(paste0(dir_move, "/trajectory_boar_9.png"))) {
walk(unique(sf_boars$id), ~plot_move_map(.x))
}
## example maps
knitr::include_graphics(paste0(dir_move, "/trajectory_boar_5.png"))knitr::include_graphics(paste0(dir_move, "/trajectory_boar_7.png"))Figure 6
dat_est_ssf <- read_rds(paste0(dir_est, "/ssf.rds")) |>
dplyr::select(id, res, fits) |>
mutate(run_id = 1:nrow(pick(id))) |>
unnest(cols = fits) |> mutate(method = "issf")
dat_est_rsf <- read_rds(paste0(dir_est, "/rsf.rds")) |>
dplyr::select(id, res, fits) |>
mutate(run_id = 1:nrow(pick(id))) |>
unnest(cols = fits) |> mutate(method = "rsf")
dat_est_wrsf <- bind_rows(
map(list.files(paste0(dir_est, "/wrsf"), full.names = TRUE), ~
read_rds(.x) |> mutate(run_id = str_extract(.x, "inst_\\d{1,2}") |> str_extract("\\d{1,2}") |> as.numeric()))
) |>
mutate(term = case_when(term == "habitat" ~ "hab",
term == "dist_water" ~ "dtw",
term == "imperviousness" ~ "imp"))
dat_est <-
bind_rows(
dat_est_wrsf |> left_join(dat_est_rsf |> dplyr::select(id, res, term, run_id)),
dat_est_rsf,
dat_est_ssf
) |>
dplyr::filter(term %in% c("hab", "imp", "dtw")) |>
mutate(
id = factor(id, levels = sort(as.numeric(as.character(unique(id))))),
sim = factor(method, levels = c("issf", "rsf", "wrsf"),
labels = c("iSSA", "RSA", "wRSA")),
c_name = factor(term,
levels = c("hab", "imp", "dtw"),
labels = c("Forest<br>habitat", "Impervi-<br>nousness",
"Distance<br>to water"))
)
## linear numeric scale with significant outcomes only
dat_est |>
filter(estimate > -20) |>
mutate(sig = fct_rev(ifelse(!(low < 0 & high > 0), "significant", "not significant"))) |>
ggplot(aes(x = res, y = estimate, group = id)) +
annotate(geom = "rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = 0, alpha = .12) +
geom_hline(aes(yintercept = 0), color = "grey67") +
geom_line(
aes(color = id, color = after_scale(desaturate(lighten(color, .33), .5))),
linewidth = .6
) +
geom_point(
fill = "white", color = "white", shape = 21, stroke = .6, size = 2.2
) +
geom_point(
aes(fill = id, alpha = sig), shape = 21, stroke = .6,
size = 2.2
) +
facet_grid(
c_name ~ sim,
scales = "free_y",
space = "free_y"
) +
coord_cartesian(clip = "off") +
scale_x_continuous(
breaks = sort(unique(dat_est$res)),
labels = sort(unique(dat_est$res)),
expand = c(.01, .01)
) +
scale_y_continuous(expand = c(.07, .07), breaks = seq(-20, 5, by = 1)) +
scale_color_manual(
values = pal, name = NULL, labels = function(x) paste0("WB", x),
guide = guide_legend(nrow = 1, order = 1)
) +
scale_fill_manual(
values = pal, name = NULL, labels = function(x) paste0("WB", x),
guide = guide_legend(nrow = 1, order = 1)
) +
scale_alpha_manual(values = c(1, .25), guide = guide_legend(override.aes = list(fill = "grey40"), order = 2)) +
labs(x = "Sampling interval in minutes", y = "Estimated selection parameter", alpha = NULL) +
theme(
axis.text.x = element_text(size = rel(.7), hjust = c(1.05, .25, .25, .4, .5, .75)),
axis.text.y = element_text(size = rel(.7)),
axis.title.y = element_text(margin = margin(r = 3)),
panel.spacing.x = unit(1.4, "lines"),
legend.position = "bottom",
legend.text = element_text(family = "PT Sans Narrow"),
legend.justification = "left",
legend.box.just = "left",
#legend.location = "plot",
legend.key.width = unit(1.4, "lines"),
legend.margin = margin(1, 1, 1, 1),
strip.text.x = element_markdown(hjust = 0),
strip.text.y = element_markdown(angle = 0, hjust = 0, vjust = .5, lineheight = .9,
margin = margin(l = 6, t = -.3)),
strip.background = element_rect(color = "transparent", fill = "transparent")
)ggsave(paste0(dir_plots, "/6_emp_est_env.png"), width = 10.5, height = 11.5, dpi = 600, bg = "white")Figure S6
Maps of the environmental variables and the study area
Plot Environmental Variables
phab <-
ggplot(df_env, aes(x, y, fill = hab, color = after_scale(fill))) +
geom_tile(size = .1) +
scale_fill_manual(
values = c("#d7b04e", "#005a00"),
name = "Forest habitat"
) +
coord_equal(expand = FALSE, clip = "off") +
guides(fill = guide_legend(title.position = "top", title.hjust = .5)) +
theme(axis.title = element_blank(), legend.text = element_text(size = 15, family = "PT Sans"))
pimp <-
ggplot(df_env, aes(x, y, fill = imp, color = after_scale(fill))) +
geom_tile(size = .1) +
ggsci::scale_fill_material(
palette = "light-green",
name = "Imperviousness (scaled)"
) +
coord_equal(expand = FALSE, clip = "off") +
guides(fill = guide_colorbar(
barwidth = unit(24, "lines"), barheight = unit(.6, "lines"),
title.position = "top", title.hjust = .5)
) +
theme(axis.title = element_blank())
pdtw <-
ggplot(df_env, aes(x, y, fill = dtw)) +
geom_tile(size = .1) +
scale_fill_carto_c(
palette = "ag_Sunset",
direction = -1,
limits = c(0, NA),
name = "Distance to water (scaled)"
) +
coord_equal(clip = "off", expand = FALSE) +
guides(fill = guide_colorbar(
barwidth = unit(24, "lines"), barheight = unit(.6, "lines"),
title.position = "top", title.hjust = .5)
) +
theme(axis.title = element_blank())
p_s6a <-
(pdtw + labs(x = "", y = "") +
pimp + labs(x = "", y = "Latitude") +
phab + labs(x = "Longitude", y = "")) *
theme(
plot.margin = margin(5, 15, 5, 0),
panel.background = element_rect(fill = "grey85", color = "grey85"),
axis.text = element_text(size = rel(.75)),
legend.title = element_text(size = rel(1.4), face = "bold")
) +
plot_layout(ncol = 1)
ggsave(paste0(dir_plots, "/_raw/s6a_maps_emp_env_vars.png"), width = 8.6, height = 18, dpi = 600, bg = "white")
p_s6aOverview Map
globe <- d6berlin::globe(size_pin = 3, col_earth = "grey96")
sf_europe <-
ne_countries(continent = "Europe", scale = 10, returnclass = "sf") |>
st_transform(crs = st_crs(sf_boars))
bbox <-
sf_boars |>
st_bbox() |>
st_as_sfc()
germany <-
ggplot(sf_europe) +
geom_sf(aes(fill = sovereignt == "Germany"), color = NA, size = 0) +
geom_sf(data = bbox, fill = "#415c27", size = 1.2, color = "black") +
coord_sf(xlim = c(3990000, 4720000), ylim = c(2660000, 3580000), expand = FALSE) +
scale_fill_manual(values = c("grey92", "#c3dbab"), guide = "none") +
theme_void()
p_s6b <- globe / plot_spacer() / germany + plot_layout(heights = c(1, .02, 1.6))
ggsave(paste0(dir_plots, "/_raw/s6b_study_area.png"), width = 9, height = 18, dpi = 600, bg = "white")
p_s6bFigure S7
Correlations of environmental variables
## extract covariates
## enough to use only 30 minutes data, because the others are subsets
df_boars_env <-
df_boars |>
filter(res == 30) |>
rename(t = "time") |>
with(track(x, y, t)) |>
extract_covariates(raster::stack(ras_env)) |>
mutate(hab = factor(hab, levels = c(0, 1), labels = c(" non-forest ", " forest")))
df_boars_env$id <- df_boars |> filter(res == 30) |> mutate(id = paste("Wild boar", id)) |> pull(id)
## compute correlation matrix
cors <-
df_boars_env |>
mutate(hab = as.numeric(hab) - 1) |>
dplyr::select(
"Forest\nhabitat" = hab,
"Imperviousness" = imp,
"Distance\nto water" = dtw
) |>
correlate(method = "spearman", diagonal = 1) |>
shave(upper = FALSE)
cors <- cors |>
pivot_longer(
cols = -term,
names_to = "colname",
values_to = "corr"
) |>
mutate(
rowname = fct_inorder(term),
colname = fct_inorder(colname),
lab = ifelse(!is.na(corr), sprintf("%1.3f", corr), "")
)
ggplot(cors, aes(rowname, fct_rev(colname), fill = corr)) +
geom_tile() +
geom_text(aes(
label = lab, color = abs(corr) < .75),
family = "PT Mono", size = 7
) +
coord_fixed(expand = FALSE) +
scale_color_manual(values = c("white", "black"), guide = "none") +
scale_fill_distiller(
palette = "PuOr", na.value = "white",
direction = 1, limits = c(-1, 1),
name = "Spearman\nCorrelation"
) +
guides(fill = guide_colorbar(barheight = unit(10, "lines"), barwidth = unit(.6, "lines"),
draw.ulim = FALSE, draw.llim = FALSE)) +
labs(x = NULL, y = NULL) +
theme(
axis.text = element_text(family = "PT Sans Narrow", size = rel(1.2)),
axis.text.y = element_text(hjust = 0),
legend.position = "inside",
legend.position.inside = c(.88, .75),
legend.title = element_text(margin = margin(1, 0, 15, 0)),
legend.text = element_text(hjust = 1)
)ggsave(paste0(dir_plots, "/s7_corr_covar_emp.png"), width = 8, height = 7, dpi = 600, bg = "white")Figure S11
Environmental variables per location
p_s16 <-
ggplot(df_boars_env, aes(dtw, imp, color = hab)) +
geom_vline(xintercept = 0, color = "grey65", linetype = "23") +
geom_hline(yintercept = 0, color = "grey65", linetype = "23") +
geom_point(shape = 21, fill = "white", stroke = 1.2) +
geom_point(color = "white") +
geom_point(alpha = .1) +
facet_wrap(~id, ncol = 4) +
scale_color_manual(values = c("#a3811e", "#228822"), name = "Forest habitat:") +
labs(x = "Distance to water (scaled)", y = "Imperviousness (scaled)") +
theme(plot.margin = margin(1, 10, 1, 1))
ggsave(paste0(dir_plots, "/s11_env_vars_move_by_id.png"), width = 15, height = 12, dpi = 600, bg = "white")
p_s16Session Info
# git2r::repository() ## remove information to ensure double-blind review
sessionInfo()R version 4.4.2 (2024-10-31)
Platform: aarch64-apple-darwin20
Running under: macOS Ventura 13.6.4
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.0
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
time zone: Europe/Berlin
tzcode source: internal
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] here_1.0.1 corrr_0.4.4 amt_0.3.0.0
[4] patchwork_1.2.0 colorspace_2.1-0 rcartocolor_2.1.1
[7] rnaturalearth_1.0.1 terra_1.7-78 sf_1.0-16
[10] purrr_1.0.2 lubridate_1.9.3 ggtext_0.1.2
[13] ggplot2_3.5.1 stringr_1.5.1 forcats_1.0.0
[16] tidyr_1.3.1 dplyr_1.1.4 readr_2.1.5
loaded via a namespace (and not attached):
[1] Rdpack_2.6 DBI_1.2.3
[3] s2_1.1.6 remotes_2.5.0
[5] rlang_1.1.4 magrittr_2.0.3
[7] e1071_1.7-14 compiler_4.4.2
[9] systemfonts_1.1.0 vctrs_0.6.5
[11] profvis_0.3.8 wk_0.9.1
[13] pkgconfig_2.0.3 fastmap_1.1.1
[15] backports_1.5.0 ellipsis_0.3.2
[17] labeling_0.4.3 utf8_1.2.4
[19] promises_1.3.0 rmarkdown_2.26
[21] markdown_1.12 sessioninfo_1.2.2
[23] tzdb_0.4.0 ragg_1.3.2
[25] rnaturalearthhires_1.0.0.9000 xfun_0.49
[27] cachem_1.0.8 jsonlite_1.8.8
[29] later_1.3.2 R6_2.5.1
[31] RColorBrewer_1.1-3 stringi_1.8.4
[33] pkgload_1.3.4 Rcpp_1.0.12
[35] knitr_1.46 usethis_2.2.3
[37] httpuv_1.6.15 Matrix_1.7-1
[39] splines_4.4.2 timechange_0.3.0
[41] tidyselect_1.2.1 rstudioapi_0.16.0
[43] yaml_2.3.8 codetools_0.2-20
[45] miniUI_0.1.1.1 pkgbuild_1.4.4
[47] lattice_0.22-6 tibble_3.2.1
[49] shiny_1.8.1.1 withr_3.0.0
[51] evaluate_0.23 survival_3.7-0
[53] units_0.8-5 proxy_0.4-27
[55] urlchecker_1.0.1 xml2_1.3.6
[57] pillar_1.9.0 KernSmooth_2.23-24
[59] checkmate_2.3.1 generics_0.1.3
[61] sp_2.1-4 rprojroot_2.0.4
[63] hms_1.1.3 commonmark_1.9.1
[65] munsell_0.5.1 scales_1.3.0
[67] xtable_1.8-4 class_7.3-22
[69] glue_1.8.0 tools_4.4.2
[71] d6berlin_1.0.0 fs_1.6.4
[73] grid_4.4.2 rbibutils_2.2.16
[75] devtools_2.4.5 raster_3.6-26
[77] cli_3.6.3 textshaping_0.4.0
[79] fansi_1.0.6 d6_0.1.0.4
[81] gtable_0.3.5 ggsci_3.1.0
[83] digest_0.6.35 classInt_0.4-10
[85] htmlwidgets_1.6.4 farver_2.1.1
[87] memoise_2.0.1 htmltools_0.5.8.1
[89] lifecycle_1.0.4 httr_1.4.7
[91] mime_0.12 gridtext_0.1.5